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When  a  great  circle  path  is  traveled  between  two  points  on  the  earth’s 
surface,  the  heading  generally  varies  in  time.  Formulas  are  presented  below 
for  the  initial  and  final  headings  on  the  shortest  great  circle  path  between  an 
origin,  O,  and  a  destination,  D,  whose  positions  are  specified  by  latitude  and 
longitude.  A  formula  is  also  given  for  the  great  circle  distance  between  O 
and  D. 

In  appendix  A,  a  Fortran  computer  routine  is  provided  for  implementing 
these  formulas  for  great  circle  headings  and  distances. 

RESULTS 

Assume  that  the  latitude  and  longitude  for  the  origin  and  destination 
are  given  as 

0  =  (LAT0,  LNGJ 
D  =  (LAT^  LNGJ  . 

It  is  assumed  that  ail  latitudes  lie  in  the  range  [ — 90°,  +  90°],  positive  being  to 
the  north.  Likewise  all  longitudes  lie  in  the  range  [  —  180°,  +180°],  positive 
being  eastwards.  One  computes  the  quantities 


ALAT  =  LATq  -  LATQ 
A LNG  =  LNGd  -  LNGQ 
LAT  =  1/2  ( LATn  +  LATn)  . 

ave  DO 

The  shortest  great  circle  distance  between  0  and  D  may  be  found  from 
the  formula 

d  =  60co8_1  IsinCLArQjsinfLAr^)  +  cos(Z.A7’0)cos(£,A7,i?)cos(ALArG)|  .  (1) 

In  equation  1,  the  arc  cosine  function  is  assumed  to  take  on  values  in  the 
range  [0°,  180°],  and  the  distance,  d,  is  calculated  in  nautical  miles. 


VVj 


Before  presenting  the  general  formulas  for  initial  heading  H0  and  final 
heading  HD,  several  special  cases  must  be  considered  that  correspond  to 
travel  along  a  meridian.  These  cases  arise  when  A LNG  =  0°  or  ±  360°, 
when  AUVG  =  ±180°,  or  when  0  or  D  is  at  the  north  or  south  pole. 

'When  O  (for  example)  is  the  north  pole,  its  longitude  is  indeterminate 
and  the  case  should  be  treated  as  having  ALNG  =  0°.  Likewise,  for  the  south 
pole  and  for  D.  If  A LNG  =  0°  or  ±  360° ,  the  rules  are: 


•  If  MAT  >  0°,  H0  =  Hd  =  0°. 


•  If  MAT  <  0°,  H0  =  Hd  =  180°. 

•  If  MAT  =  0°,  no  movement  is  required,  and  the  headings  are 
indeterminate. 

On  the  other  hand,  if  A LNG  =  ±  180°,  travel  over  the  north  or  south 
pole  is  generally  required.  Therefore,  the  rules  are: 

•  If  LATau  >  0°,  the  path  crosses  the  north  pole,  so  that  HQ  =  0° 
and  Hd  =  180°. 

•  If  LATave  <  0°,  the  path  crosses  the  south  pole,  so  that  HQ  -  180° 
and  Hd  =  0°. 

When  LATaut  =  0°,  the  headings  are  indeterminate,  since  O  and  D  are  on 
diametrically  opposite  sides  of  the  earth.  For  any  initial  heading  HQ,  a  great 
circle  path  may  be  found  to  D,  and  it  will  have  HD  ~  HQ. 

When  the  foregoing  special  cases  do  not  apply,  the  initial  heading  may 
be  found  from  the  formula 


cos  (AL  AT/2) 

tan  (ALNGI2)  sin  (LAT  ) 

aue 


-  Tan~l 


sin  (AL  AT/2) 

tan  (ALNGI2)  cos  (LAT  ) 

aue 


+  180°-  180° 


sign  (tan  (ALNGI2))  if  LAT ^  >  0° 
0  if  LAT  <5  0° 

out 


(2) 


In  equation  2,  Tan~l  represents  the  principal  determination  of  the  arc 
tangent  function,  taking  on  values  in  the  interval  (-90°,  90°).  The  final 


heading  is  given  by 


Hd  =  -Tan~l 


cos  (AL  AT/2) 

tan  (ALNG/2)  sin  (LAT  ) 

ave 


-Tan'1 


sin  (AL AT/2)  . 

tan(ALNG/2)  cos  (LAT  ) 

ave 


+  180°  -  180° 


0  if  LAT  a  0° 

ave 

sign  (tan  (ALNG/2))  if  LAT ^  S  0° 


(3) 


DERIVATION  OF  HEADING  FORMULAS 

The  geometry  used  to  derive  equations  1  through  3  is  shown  in  figure  1. 
Two  cases  are  illustrated:  an  "eastbound”  case,  in  which  A LNG  z  (0°,  180°)  or 
equivalently  A LNG  e  (—360°,  —180°),  and  a  "westbound”  case,  in  which 
A  LNG  z  (  —  180°,  0°)  or  A LNG  e  (180°,  360°).  In  both  cases,  the  path  to  be 
traversed  is  marked  with  arrows.  Meridians  have  been  indicated  from  the 
north  pole  to  O  and  D.  The  lengths  of  these  segments  of  meridians,  a  and  6, 
and  the  length,  c,  of  the  great  circle  path  will  be  measured  as  angles 
subtended  by  the  arcs  at  the  center  of  the  earth.  The  internal  angles  A,  B, 
and  C  of  the  spherical  triangles  have  been  indicated. 

In  Case  I,  the  eastbound  case,  we  have 


Hd  =  180°  -  B  ,  (4) 

while  in  Case  II,  the  westbound  case, 


The  angle  C,  which  is  measured  at  the  north  pole,  may  be  written  in  the  form 


_  |  |ALM?|  if  |ALM?|  <180° 

C  l  360°  -  |ALM?|  if  |ALNG|  >  180°  *  (6) 

so  that  C  is  always  positive  and  less  than  180°. 

Equation  1  for  the  great  circle  distance  is  easily  obtained  from  the  law  of 
cosines  for  spherical  triangles:1 


cos(c)  =  cos  (a)  cos  (6)  +  sin  (a)  sin  (6)  cos(C)  . 


When  a  and  b  are  expressed  in  terms  of  LATD  and  LAT0,  and  equation  6  is 
used,  the  result  is 

c  =  cos-1  ^sin(LATQ)  sin(LAT ^  +  cos(LATQ)  cosiLAT  Q )  cos(ALM?)|  . 

One  degree  of  arc  subtended  at  the  earth’s  center  corresponds  to  a  distance  of 
60  n.mi.  on  the  earth’s  surface,  so  that  equation  1  immediately  follows. 

It  is  clear  from  equations  4  and  5  that  the  desired  headings  may  be  found 
once  A  and  B  are  known.  The  spherical  triangles  in  figure  1  may  be  solved 
to  give  the  angles  A  and  B  in  terms  of  a,  6,  and  C,  using  Napier’s 
analogies  and  standard  rules  for  determining  quadrants.  The  appropriate 
formulas  are 


sin  (a  —  6)  /  2  _  tan  (A  —  B)  12 
sin(a  +  b)/ 2  cotC/2 


cos  (a  —  b )/ 2  _  tan  {A  +  B)/ 2 
cos  (a  +  b)  12  cot  Cl  2 


(8) 


The  quantities  a  —  b  and  a  +  b  may  easily  be  related  to  LATavt  and  A LAT. 
Furthermore,  equation  6  relates  C  to  | ALNG\.  Thus,  equations  7  and  8  may 


1.  Taken  from  Chemical  Rubber  Co.,  CRC  Standard  Mathematical  Tables,  20th  edition, 
1972,  pp.  200-201. 
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A  -  B 


=  ton 


-l 


—sin(ALAT/2) 


i, 

*7 


A  +  B 


=  ton-1 


torejAZ<#G/2|  sin  (LAT  ) 


cos(ALAT/2) 
ton|ALA/G/2|  sinCLAT^) 


(9) 


(10) 


Note  that  the  inverse  tangent  functions  appearing  in  equations  9  and  10  are 
not  necessarily  the  principal  determination  of  arc  tangent.  Another  way  of 
stating  this  is  to  say  that  (A  —  B)/2  and  (A  +  B)/2  are  determined  by  these 
equations,  apart  from  possible  addition  or  subtraction  of  180°. 

The  task  now  is  to  resolve  the  quadrants  in  which  the  angles  lie.  The 
standard  rules  for  doing  so  are  stated  in  CRC  Standard  Mathematical  Tables 
as: 


RULES  FOR  DETERMINING  THE  QUADRANT 
OF  A  CALCULATED  PART  OF  AN  OBLIQUE 
SPHERICAL  TRIANGLE 

(a)  If  A  >  B  >  C,  then  a  >  b  >  c. 

(b)  A  side  (angle)  which  differs  more  from  90°  than  does  another  side 
(angle)  is  in  the  same  quadrant  as  its  opposite  angle  (side). 

(c)  Half  the  sum  of  any  two  sides  and  half  the  sum  of  the  opposite 
angles  are  in  the  same  quadrant. 


Rule  (c)  implies  that  (A  +  B)/2  is  in  the  same  quadrant  as  (a  +  b)/2 
90°  —  LAT  Thus,  one  has 


0°  <  1/2  (A  +  B)  s  90°  Lf  LAT  2  0° 

ave 


90°  <  1/2  (A  +  B)  <  180°  if  LAT  <  0° 

ave 
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Consequently,  equation  10  gives 


A  +  B 


=  Tan' 


cos (ALATI2)  1  [  °°  ifLAroU*a0° 

tan\ALNGI2\(LAT  )]  +  l  180°  if LAT  <0°  ’ 


where,  as  noted  earlier,  the  principal  determination  of  arc  tangent,  Tan~l, 
takes  on  values  between  —90°  and  +90°. 

Notice  that  by  construction,  A  and  B  both  lie  in  the  interval  (0°,  180°). 
(The  limiting  cases  in  which  A  or  B  equals  0°  or  180°  have  been  handled 
separately.)  It  follows  that 

(A  -  B)/2  e(— 90°,  90°)  . 

As  a  result,  the  inverse  tangent  function  in  equation  9  must  be  identified  as 
Tart-1.  Equations  9  and  11  may  then  be  combined  to  obtain  expressions  for  A 
and  B.  When  these  are  inserted  into  equations  4  and  5,  cases  I  and  II  may  be 
combined  into  the  expressions  for  HQ  and  HD  in  equations  2  and  3. 


that: 


By  examining  various  limits  of  equation  2,  it  is  straightforward  to  check 


•  H0  is  continuous  at  ALNG  =  0°,  for  ALAT  *  0°. 

•  H0  is  continuous  at  ALAT  —  0°,  for  ALNG  *  0°. 

The  heading  may  take  on  any  value  between  0°  and  360°  if  ALAT  and  ALNG 
are  allowed  to  go  to  zero  simultaneously.  This  is  the  mathematical  realization 
of  the  earlier  observation  that  heading  is  indeterminate  in  this  limit. 

Similarly,  one  finds  that: 


•  H0  is  continuous  at  LATavt  =  0,  for  ALNG  *  ±180°. 

•  H0  is  continuous  at  ALNG  ~  ±  180°,  for  LATaue  *  0°. 

Again,  a  singularity  is  present  at  the  point  LATave  -  0°,  ALNG  =  ±  180°. 
Any  value  of  H0  may  be  obtained  as  this  point  is  approached,  depending  on 
the  limiting  relation  between  LATaue  and  ALNG.  This  reflects  an 
indeterminacy  in  heading  mentioned  above. 
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APPENDIX  A 


LISTING  OF  THE  FORTRAN  SUBROUTINE  GRCRCL 


SUBROUTINE  GRCRCLC  LATC.  LONG »  LATO,  LQNO  »  OIST,  HSAONG  ) 

CCMPUTS  THE  01  STANCE  (OIST  CMM>)  ANO  INITIAL  HEADING  (H  =  A0NG 
COEG))  BRGM  POSITION  LAT0-L0N0  COECIMAL  OEGRSES)  TQ  POSITION 
LATO-LONO. 

SEAL  LATO  .LONO, LATO, LQNO, 01 ST ,HE AONG 

REAL *9  PI, OcGRAD.EPS, HALF 

PARAMETER  CPI  *  3*1415 9 265353979300) 

PARAMETER  COEGRAO  «  PI/130.000) 

PARAMETER  CEPS  »  t.00-7) 

PARAMETER  CHALF  a  5.03-1) 

REALMS  OLATC.OLONO,  CL  ATO  ,  OLQNO  ,  HOLL  ON,  HOLLA  T 

5£AL*9  A VGL  AT , RAOHOG , TEMPI , TEMP  2 

REAL**  0 A3S, OS  IN, OCOS.OT AN, OAS I N  ,0AC0S , 0 ATAN 

OLATG  a  0£GRA0a09LsCLAT0) 

OLONO  a  0  EGR A04D3LSCL2N3) 

0LAT3  a  0EGRA0*03L=CLAT0) 

OLONO  a  0EGRA0*03LECL0N0) 

HCLLQN  a  C OLQNO  -  0L3H0)*HALF 
10  IF  C  HOLLO  N  .LS.  -PI*HALF)  THEN 
HOLLON  a  HOLLON  ♦  PI3HALP 
GO  TO  10 
ENOIF 

20  IF  CHOLLON  .GT.  ?I*MALF)  THEN 
HOLLON  a  HOLLON  -  PI*mALP 
GO  TO  20 
ENOIF 

-tOLLAT  a  COLATO  -  3L  AT  G)  PHAL  F 
A VGL  AT  a  COLATO  ♦  3L  ATO ) PHALF 


OIST  a  60 .0O0*OACCSC3SlNCOLATO)*OSINCOLATO)>CCGSCOLAT3)* 
I  0CO3C0LAT0)4OCGSCQL3MQ-QLQNQ))/0EGRAQ 

SPECIAL  CASES. 


I 

I 


IP  CC:A3SCH0LL0N).LS.SPS).3R .COAaSCHOLLQN).GE.PI*HALF-SPS)) 
I  THEN 

H E  AONG  «  P I/O SGR  AO 

Is  CAVGLAT.GE. 0.000)  HEAOnG  »  0.0E0 
RETURN 
ENOIF 

IF  C  OAES (rtOLLAT).Gc.P C*hALb-EOS )  THEN 
OIST  a  0.3E3 
1; AONG  a  0.0 eO 
RETURN 
ENOIF 


ooo 


SEGULA3  CASE. 

TEMPI  «  3aTiM<3C3S(H0U.AT)/C0TANCH0U.0N)*0SIN<AVGLAT)>) 
TE.MP2  *  .  0ArAN(3S:M<H3LUAn/C0T»*l<H0LUaM)*0C3SCAVGLAT)>) 
5A0H3G  *  T;MP1-T i-P ’ *PI 

IP  CAVGLAT  .Ge.  0.000)  PAOHOG  *  SA0H0G  -  0  Si GnC P I ,riG U.QN) 
MEAOMG  *  OA0M3G/05G? AO 


setup* 

cN3 


A- 2 


